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We investigate cage breaking in dense hard disk systems using a model of three Brownian disks 
confined within a circular corral. This system has a six-dimensional configuration space, but can be 
cquivalently thought to explore a symmetric one-dimensional free energy landscape containing two 
energy minima separated by an energy barrier. The exact free energy landscape can be calculated 
as a function of system size. Results of simulations show the average time between cage breaking 
events follows an Arrhenius scaling when the energy barrier is large. We also discuss some of 
the consequences of using a one-dimensional representation to understand dynamics in a multi- 
dimensional space, such as diffusion acquiring spatial dependence and discontinuities in spatial 
derivatives of free energy. 

PACS numbers: 64.70.pm,65.40.gd,71.55. Jv 



I. INTRODUCTION 

The concept of a potential or free energy landscape 
for condensed matter systems is an appealing one [IH3], 
and its use to understand the glass transitions in vari- 
ous materials has received considerable attention |3Ml3|. 
Each configuration of a system can be assigned an en- 
ergy based on how individual elements interact, and the 
set of all configurational energies then constitutes an en- 
ergy surface or "landscape" . The temporal evolution of 
a system can be thought of as an exploration of various 
topographical features on the energy surface, with re- 
laxation events understood as motions between adjacent 
energy minima [5J rBHT!)] . As the glass transition is ap- 
proached, the system finds itself in progressively deeper 
local minima on the landscape until thermal energy is un- 
able to excite the system over a barrier and into a lower 
energy state. 

A complete energy landscape accounts for all degrees 
of freedom within a system, each of which contributes a 
single dimension to the configuration space. For many 
body systems then, the complete energy landscape can 
be extraordinarily complex. A central problem in using 
the energy landscape approach is minimizing the num- 
ber of "reaction coordinates" (likewise "order parame- 
ters") while still adequately describing the surface [TBl - 
[23] . Studies often choose one or two coordinates of inter- 
est and examine how the system evolves on representative 
ID or 2D landscapes, though in some cases such low di- 
mensional projections may introduce artifacts or occlude 
important features of the landscape |2"TH2"3"] . 

These complexities aside, the typical picture of an en- 
ergy landscape is one of a surface containing hills and 
valleys. However, for the case of a purely hard core po- 
tential, such as that of hard disks or hard spheres, this 
picture is incorrect. In a hard disk or hard sphere sys- 
tem, the Hclmholtz free energy is governed entirely by 
entropy, 

F = -TS, (1) 

to within an additive constant (25] [26] . As all allowable 
configurations have identically zero potential energy, they 




FIG. 1. (Color online) (a). Three distinguishable hard 
disks are confined within a circular corral of fixed size. The 
variables h and 6 are macrostate variables, to be defined 
later. Dashed lines are lines of constant h and 0. (b). A 
2D slice through the 6D configuration space with constant 
(rL,r2). The yellow (gray) region is accessible to the center 
of disk 3, whereas the dark gray region is energetically forbid- 
den. Upper and lower dashed lines represent possible system 
macrostates at fixed h, where the of the line length determines 
the entropy of the state. 

are all equiprobable. Thus, in terms of the total config- 
uration space, the energy landscape is completely flat. 
Upon projection into a lower dimensional space, however, 
hills and valleys related to entropic minima and maxima 
can arise [5J E H EZH25]. 

II. MODEL SYSTEM 

Consider a system of N = 3 hard disks confined within 
a hard circular corral, as illustrated in Fig. [lja). The 
configuration space of this system is six-dimensional (6D) 
and is completely described by all allowable combinations 
of (r\,r2,rs). However, as stated previously, the energy 
landscape is flat because each configuration is as proba- 
ble as any other. Now, we consider a 2D slice though the 
landscape examining all configurations of r*3 while (fi, ?%) 
are held at constant values, such as shown in Fig. [ijb). 
All allowed configurations are again equiprobable and are 
indicated by the light gray region (yellow online). The 
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landscape in this scenario remains flat. However, there 
are more configurations where disk 3 exists along the up- 
per dashed line than on the line between disks 1 and 2. 
If we interpret these two lines as macrostates of the sys- 
tem, the entropy of the upper line is larger and therefore 
the free energy in Eqn. ([I]) is lower. By considering all 
macrostates, that is all lines parallel to those in Fig. |ljb) , 
we obtain a ID representative landscape with variations 
in free energy. 

We can immediately relate these ideas to a dynamic 
version of the system in Fig.[lja), where all disks wander 
around the corral by means of Brownian motion. Here, 
we may ignore degrees of freedom associated with mo- 
mentum and focus only on the six spatial components. 
As the system explores 6D configuration space, disk 3 
spends more time in configurations with high entropy on 
the ID landscape, such as along the upper dashed line 
in Fig. |ljb) . To transition into the lower cusped region, 
it must pass through an entropic bottleneck. Of course, 
this behavior is not unique to disk 3 but applies to all 
disks in our system. Hence, the transition of any sin- 
gle disk corresponds to the entire system crossing a free 
energy barrier. 

Prior to a barrier crossing event, the motion of each 
disk is localized by the presence of other disks and the 
wall. As with individual molecules in dense liquids 
or particles in a dense colloidal suspension, the disks 
can then be described as "caged" [H I5UP0]. Devi- 
ations from strongly localized behavior, such as those 
during a crossing event, are considered "cage break- 
ing" (37] [33 Hi 46 . One can consider both caging and 
cage breaking from the point of view of an energy land- 
scape, where caging is the motion of the system around 
some local minimum and cage breaking is the relaxation 
over an energy barrier [29, 47 52 . 

In the spirit of [27l EH H2 EH [53] , we introduce the 
minimal model system illustrated in Fig. [l^a) in order 
to explore the relationship between caging and the free 
energy landscape of a system with hard core potentials. 
Here, a cage breaking event can be simply described as 
one disk passing between the other two. One can imagine 
for a very large corral, the motions of any disk would only 
rarely be influenced by the other disks or the boundary, 
and so dynamics would be similar to those in a dilute 
colloidal suspension. However, this picture changes for 
smaller corral sizes, or similarly higher packing fractions, 
where interactions between the disks or with the wall 
are more frequent and dynamically restrictive. Though 
somewhat contrived, these purely geometric constraints 
are of similar character to those encountered by real par- 
ticles in a densely packed system, such as colloidal super- 
cooled liquids or glasses [51H5T] . 
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FIG. 2. Trajectories through /i-space for (top-bottom): e = 
0.10, 0.20, and 0.61. Dotted lines pass through h = and 
dashed lines pass through h — ±2. The length of time be- 
tween crossing events is clearly a sensitive function of e. Small 
fluctuations about h — (such as those in the boxed regions 
for e = 0.20) are not considered as true cage breaking events. 



Choosing r = 1, we note that all lengths are by defini- 
tion in units of a particle radius. The minimum system 
size that permits cage breaking occurs at e = 0, where 
the three particles exactly fit across the diameter of the 
corral. During each simulational run, we fix the value 
of e and allow particles to execute Brownian motion, de- 
scribed below. Each run consists of 10 8 Monte Carlo 
steps (mcs) for a particular value of e, during which each 
particle has the opportunity to make a displacement. 
Values of e are chosen within a range of e G [0.037, 10.0]. 

To simulate Brownian dynamics, displacements for 
each particle at every time step in the x— and 
y— directions are sampled from a Gaussian distribution 
with variance a 2 — 2D = 1 x 10 -3 . A displacement 
is accepted if it does not result in particle-particle or 
particle-wall overlap, otherwise the offending particle re- 
mains fixed for that time step. This results in an RMS 
displacement of approximately 2% of a particle diam- 
eter at each time step. The value of D, and therefore 
the temperature, is constant across all simulations. The 
order in which particles are sampled is randomized at 
each time step and satisfies detailed balance [58]. For 
the stated value of D, the fraction of accepted displace- 
ments range from 0.954 at e = 0.037 to 0.997 at e = 10.0. 



IV. DYNAMICS IN ID 



III. SIMULATION DETAILS 

Three hard-disks of radius r = d/2 = 1 are confined 
to a hard-walled, circular corral of radius Rc = 3 + e. 



To study the dynamics of this system, we first sim- 
plify the 6D configuration space by projecting down to 
one dimension. We define a macrostate variable h, shown 
in Fig. [T|a), as the distance of the center of disk 3 nor- 
mal to a line segment drawn from the center of disk 1 
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FIG. 3. (Color online) Mean-square displacement in /t-space 
for systems with (bottom to top, dark gray to light gray) e = 
0.045, 0.06, 0.10, 0.15, 0.25, 0.40, 0.61, 1.0, 2.0, and 4.0. 
Dashed line has a slope of one. 



to that of disk 2. Put another way, if a line drawn from 
the center of disk 1 to that of disk 2 defines the pos- 
itive IE-axis of a coordinate system, h is given by the 
y-coordinate of disk 3. Therefore, h can take positive or 
negative values, which from geometry are limited to the 
range [-h max , h max ], where h ma , x = 2 + e + V'S + 4e + e 2 . 
This definition allows for a cage rearrangement to be de- 
scribed as the system passing through h = 0, regardless 
of which particle passes between the others. Finally, we 
note that this definition maps all rotationally symmetric 
states to the same value of h. 

Shown in Fig.[2]are trajectories in ft,-space for three val- 
ues of e. As e increases from top to bottom, the length 
of time between cage breaking events decreases, as one 
would expect. For these and all other trajectories, we 
find systems spend the majority of time localized around 
h = ±2, commensurate with a particle diameter. The 
strength of the localization can be inferred from the fluc- 
tuations about h = ±2, which increase with increasing e. 
Hence, as the system size increases, the strength of the 
localization decreases and cage breaking events become 
more frequent. 

To quantify motions in /i-space, we define a 
one-dimensional mean-square displacement (MSD), 
(Ah 2 (At)) = ([h(t + At) - h(t)f) , where the average 
extends over equivalent lag times At. Shown in Fig. [3] 
are (Ah 2 (At)) for various e. For all system sizes, we 
observe diffusive behavior on short time scales such that 
(Ah 2 (At)) oc At, however effects of finite system size 
are apparent from differences in the intercept at At = 1, 
as will be discussed shortly. 

For the largest e, (Ah 2 (At)) retains diffusive behav- 
ior until finally plateauing due to the finite size of the 
system. However, the onset of distinctly subdiffusive be- 
havior and development of increasingly long plateaus are 
apparent as e is decreased. For the smallest values of e in 
Fig. [3j (Ah 2 (At)) is qualitatively similar to the MSDs of 



supercooled liquids and glass formers, where the plateau 
indicates timescales over which caging occurs. The up- 
turn of (Ah 2 (At)) indicates that cage breaking eventu- 
ally occurs for all systems. 
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FIG. 4. (a) Short time diffusion coefficient scaled by squared 
particle size, (b) Short time diffusion coefficient scaled by 
squared system size, (c) Average time between cage breaking 
events. Error bars are calculated from the statistical uncer- 
tainty based on the number of events. Where not visible, 
error bars are smaller than the symbol, (d) Height of the 
free energy barrier as a function of e. Circles are the results 
from simulations and the solid line is from the calculations 
described in Sec.jv] Dashed line grows as In [e~ 7 ^ 2 ] and shows 
scaling behavior as e — > 0. 

Shown in Fig. [4^a,) are the short time diffusion coeffi- 
cients in /i-space, defined as, 



D h = 



(Ah 2 (At = l)) 



(2) 



In Cartesian space, the short time diffusion coefficient is 
a constant with no spatial or system size dependence. In 
our ID coordinate system, however, the diffusion coeffi- 
cient becomes a function of space (discussed in Sec. M) 
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and system size. For small systems, Dh is relatively con- 
stant, though we do find a weak linear dependence on 
/i max (not shown). As e increases, however, Dh increases 
without bound, again only as a consequence of our defi- 
nition of h. In large systems, two particles may be near 
each other, while the other is some distance away. In this 
case, a small displacement of one of the proximal parti- 
cles in Cartesian space can translate into a very large 
displacement in /i-space. Hence, the unbounded increase 
of Dh with e is not surprising. If the diffusion coefficient 
is normalized by the squared system size, however, as in 
Fig. |4^b), we find a trend opposite to that seen in (a). 
These values at small e are again relatively constant, but 
decrease significantly as the system size increases. Thus, 
as e increases, the system explores absolute space more 
quickly, but relative to the size of the system, the explo- 
ration of space occurs slower. 

To quantify a relaxation time scale, or the time the 
system spends caged, we define a transition time scale r 
as the average time needed to cross h = from positive 
h to negative h, or vice versa. As highlighted in Fig. [2] 
there can be small fluctuations about h = that are not 
true cage breaking events. To minimize biasing r toward 
lower timescales due to this sort of rattling, we stipulate 
that once h = is crossed, the system must move a fur- 
ther distance I* before returning. Otherwise, no crossing 
event is registered. If F(h) defines the free energy land- 
scape (see Sec. [V] for more details), the distance I* is 
calculated as 



^ J h 2 F{h) dh 

-2 

/ F{h) dh 

\ -2 



1/2 



(3) 



that is, the standard deviation of the energy landscape 
between [—2, 2]. As will be discussed in Sec.|Vj the height 
of the energy barrier decreases with system size but ex- 
tends over exactly this domain for all systems. Hence, 
Eqn. (|3) provides a consistent length scale for where the 
energy barrier drops to « 3/5 of its peak value (further 
details of calculating F B are given in Sec. [V]). Shown in 
Fig. |4^c) , as e decreases from 1 toward 0.1, the relaxation 
time increases, by approximately two orders of magni- 
tude, and continues to grow dramatically as the system 
becomes smaller. The smallest relaxation time occurs 
when e w 3.0. For e > 3.0, motion through the landscape 
is limited only by diffusion and so the relaxation time 
increases as a result of the large system size. 



V. ENERGY LANDSCAPE 

From Eqn. ([I]), determining the energy landscape 
is only a matter of calculating the entropy for each 
macrostate h. For simplicity, we set ksT = 1 and write 




FIG. 5. (Color online) Comparison between theoretical en- 
ergy landscapes (solid lines) and results from simulations 



(symbols). Circles: 
e = 1.00. 



0.06, squares: e = 0.25, triangles: 



the free energy of a state relative to the ground state as, 

F(h) = -T[S(h)-S ] 

= -]n[n(h)/no], (4) 

where subscripts refer to the ground state at a given e, 
and n(h) is the number of states that map to the same h. 
We calculate the number of states or multiplicity n(h) by 
integrating over the space of allowed configurations f2 of 
the three disks while maintaining a fixed h. In general, 
this can be written 



n(h) = J dr x dr 2 dr 3 S[h - H{f\, f 2 , r 3 )]. 



(5) 



where the function H(ri,r2, 7*3) calculates the value of 
h given the coordinates of the three disks, and S is the 
Dirac delta function. The expression in Eqn. ([5]) can be 
reduced to a ID integral which we integrate numerically. 
We then calculate F(h) using Eqn. Q. Further details 
of calculating n(h) are given in the appendix. 

To determine the distribution n(h) from simulations, 
we directly count the number of states by constructing 
histograms of h with small bins [59]. At this point, F(h) 
can be computed directly using Eqn @. 

In Fig. [5] we compare analytical calculations of energy 
landscapes (solid lines) to those resulting from simula- 
tions (symbols). In all cases, the agreement between the- 
ory and simulation is excellent and so we are confident 
that the number of simulation steps adequately samples 
the configuration space. Energy landscapes for all values 
of e are symmetric about h = and are double- welled, 
with infinitely high barriers at h = h max , corresponding 
to particles being unable to escape the corral. The height 
of the energy barrier Fb, shown in Fig. (4^d) , increases as 
e decreases and diverges as e — » 0. As e — > 0, the height 
of the barrier grows as In [e~ 7 ^ 2 ] , which can be predicted 
analytically and is described in the appendix. 
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FIG. 6. Symbols: Measured relaxation time as a func- 
tion of theoretical free energy barrier height. Solid line is 
378 exp {PF B ), with /3 = 0.975 ± .04, and comes from a 
weighted fit to the open symbols. To within statistical er- 
ror, transition times grow Arrheniusly with the height of the 
energy barrier when Fb > 7.0. 



The most probable h for any corral size occurs at h = 
±2, as indicated by locations of the minima in Fig. [5] 
Inspecting Fig. [ijb), we see that there are more ways 
to place disk 3 at a value of h corresponding the cusp 
of the allowed region than anywhere else. For \h\ < 2, 
the number of states is limited by configurations where 
particles overlap, whereas this constraint disappears for 
\h\ > 2. 

Given the energy landscape, we can now ascribe the 
localization about h — ±2 seen in Fig. [2] to the system 
being trapped in one of two local energy minima. Fur- 
thermore, we see in Fig. [5] that the landscape broadens 
outward as e increases, which explains the increase in 
fluctuations shown in Fig. [2] Inspecting Fig. [3j we see 
that a longer plateau in (Ah 2 ) equates to a larger free 
energy barrier, and thus corresponds to the system being 
constrained at low e and short times to explore only those 
regions near a minimum in the energy landscape. Given 
enough time, a large thermal fluctuation allows the sys- 
tem to cross the energy barrier, producing the upturn in 
(Ah 2 (At)). 

We are also now in a position to relate the energy 
landscape to the previously measured relaxation times. 
Shown in Fig. |6j we find that these relaxation times scale 
Arrheniusly with barrier height, to within statistical er- 
ror, when Fb > 7.0. As shown in Fig. [2j fluctuations 
about \h\ = 2 grow with increasing e. Additionally, when 
e is small, the energy barrier is large. Hence, for large 
values of F B , the system finds it difficult to explore large 
values of \h\ and the dominating factor in relaxation is 
the height of the energy barrier. 

In contrast for larger e, the energy barrier is smaller 
and fluctuations about \h\ —2 are much more signifi- 
cant. In these cases, the system finds it easier to wander 
to larger values of \h\ and must first diffuse toward the 
barrier before crossing. This accounts for the deviation 
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FIG. 7. (Color online) First derivatives of F(h) for e = 0.10 
(solid line), 0.31 (dotted line), and 0.62 (dashed line). Dis- 
continuities exist at h = ±2. Additionally, a discontinuity 
exists in the second derivative at h = ±v3, as indicated by 
the arrow. 

from Arrhenius behavior toward longer relaxation times 
at smaller Fb- When Fb < 1.0, relaxation becomes es- 
sentially independent of the barrier height and instead 
depends only on diffusion and system size. 

In the overdamped limit of Brownian motion, Kramers 
rate theory [50] states that the time scale for crossing an 
energy barrier behaves as [STJ 

r cx (DuJmUJb)- 1 exp (AE/k B T), (6) 

where D is the diffusion constant, u) m and w& are the 
curvature at the minimum and barrier, respectively, and 
AE 3> ksT is the height of the energy barrier. In our 
case, to m is not defined. The energy landscapes shown 
in Fig. [5] are everywhere continuous, however they are 
not everywhere differentiable. In Fig. [7J we show the 
first derivative of F(h) for e = 0.10, 0.31, and 0.62. As 
a consequence of confinement and the projection in h- 
space, a kink in the free energy curve arises for all systems 
at h = ±2, and thus the first spatial derivative of free 
energy is discontinuous. The origins of the discontinuities 
are discussed in the appendix. Though the curves are 
steep, we point out that all values of dF/dh in Fig. [7] are 
finite. The discontinuity in dF/dh shows that uj m is not 
defined and so the minimum cannot be approximated 
as harmonic. Thus, simple expressions from Kramers 
rate theory are not able to predict the transition times 
shown in Fig. [6j Additionally for all systems, we observe 
a kink in dF/dh at h = ±Voj and therefore the second 
derivative of the free energy d 2 F/dh 2 is discontinuous at 
these points. 

As mentioned previously, the diffusion coefficient in 
/i-space is spatially dependent. We calculate a spatial 
diffusion coefficient, Dh(h) as in Eqn. |2) over a range 
[h — Sh, h + Sh], where Sh is /i max /100. Near h = ±2, 
we measure the limiting value of the variance up to 
but not including these points. Shown in Fig. [8] for 




FIG. 8. (Color online) Spatially dependent diffusion coeffi- 
cients for (bottom to top) e = 0.31, 0.70, and 1.40 for h > 0. 
Inset: a closer view of region near the kink in F(h). 



e = 0.31, 0.70, and 1.40, the spatially dependent dif- 
fusion coefficients appear to be continuous over all per- 
mitted values of h. To within our resolution, however, 
we cannot discern whether or not Dh(h) is everywhere 
differentiable. 

Our focus thus far has used h as a means of capturing 
the relevant cage breaking dynamics, but of course, one 
can imagine many options when reducing the configura- 
tion space from six dimensions down to one. For consis- 
tency, it is worth checking that the results shown above 
are not dependent on the choice of variable. We thus de- 
fine a second variable 9, shown in Fig. [T] that describes 
the angle between vectors r^i = r% — r\ and r^x = fa — f\. 
A cage breaking event can then be defined when the sys- 
tem crosses 9 = or 9 = ±tt. Free energy landscapes 
in #-space are shown in Fig. ^a.) for e = 0.10, 0.31, and 
0.62 [53]. We find the energy landscape is again double- 
welled, reflecting the two possible caged states, and has 
even symmetry about 9 = between [— 7T,7r]. Here, un- 
like F(h), we find that the locations of the minima in the 
energy landscape are not constant but move to smaller- 
values of 9 as e increases. For example in Fig.|9ja), min- 
ima shift from 9 = ±.245?r at e = 0.10 to 9 = ±0.225?r 
at e = 0.62. 

As shown in Fig. [9]ja), there are two energy barriers in 
#-space as opposed to the the single barrier in /i-space. 
While the representation has changed from variable h to 
9, the fundamental problem of caging has not, i.e. cage 
breaking is equally difficult or equally easy independent 
of the representation. Hence, for these measurements 
of energy to be consistent, there must be a relationship 
between the energy barriers measured in the two different 
coordinate systems. 

In Fig.^b), we show the height of the free energy bar- 
riers 0-space e as a function of the previously measured 
barrier heights in /i-space. The heights of the barriers 
are determined from the same simulation data that was 
used to calculate F(h). To calculate Fg(9), we consider 
the average probability that the system is poised to cage 



9/n 



(b) 



12 
10 

8 
6 
4 
2 



10 12 



F B (h) 



FIG. 9. (Color online) (a) Free energy landscape in terms of 9 
for e = 0.10 (solid line), 0.31 (dotted line), and 0.62 (dashed 
line), (b) Heights of energy barriers for all e in S-space versus 
the barrier height in /i-space, calculated as described in the 
text. Error bars indicate the uncertainty in local quadratic 
fits of n(0) and no. Dotted line is Fb(0) = Fs(h). 



break, either at 9 = or 9 



This yields 



F B (9) = -log{[n(0) + n(ir)}/{2n }}. 



(7) 



As shown in the appendix, n(0) = 4n(7r). Using this 
fact, we only measure n(0) and no to compute Fb{9) in 
Eqn. 0. The heights of the energy barriers in (9-space 
shown in Fig. [9jb) are in excellent agreement with those 
using h as a coordinate. The consistency between these 
two measurements demonstrate that both coodinates ad- 
equately describe caging and that the measured free en- 
ergies are indeed those relevant to cage breaking in our 
system. 

As in the case of F(h), F(9) is everywhere continu- 
ous. We show in Fig. 10 the first derivative of F(0). 
Here, curves are obtained by numerically differentiating 
those in Fig. [9]Ja). In contrast to F(h), dF/d9 is also 
everywhere continuous, and so F(9) is also everywhere 
differentiable. Without an analytical calculation of F(9) 
we cannot draw definitive conclusions about the conti- 
nuity of d 2 F/d9 2 . However, we do find for each e the 
hint of a subtle kink in dF/d9 located at 9 = ±7r/3, as 
well as a second e-dependent region where the slope of 
dF/d9 changes sharply, e.g. at approximately 9 = 0.167T 
for e = 0.10 in Fig. [TUj As in /i-space, diffusion in 
0-space also exhibits a complicated spatial dependence 
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FIG. 10. (Color online) First derivative of F(0) for e = 0.10 
(solid line), 0.31 (dotted line), and 0.62 (dashed line). All 
first derivatives exhibit a subtle kink at 8 — ±7r/3 indicated 
by the arrow. Horizontal dashed-dotted line is at dF/d9 = 0, 
and helps illustrate the gradual movement of the minimum to 
smaller values of 9 for increasing e. 



(not shown). One assumption leading to Eqn. ([6| is that 
D is spatially independent. Predictions for r which in- 
corporate spatially dependent diffusion are highly non- 
trivial [52] . 



successfully in /i-space or (9-space demonstrates that pro- 
jections of the landscape to a lower dimension are not 
necessarily unique. 

Depending on the choice of coordinates, the energy 
landscape of hard disk and hard sphere systems may be 
locally non-differentiable. This observation highlights the 
notion that projection of a many dimensional landscape 
down to a lower dimensional representation may intro- 
duce artifacts [2"TH2"3"] . For systems with N > 3, it is 
probable that discontinuities in spatial derivatives of free 
energy disappear completely, or perhaps are only present 
for higher order derivatives. 

For any purely soft core potential, the observed dis- 
continuities will disappear. For a potential more repre- 
sentative of colloidal suspensions, such as a hard core 
potential with short ranged repulsion |64L 165] , we expect 
the essence of our results to be valid, but in need of some 
modification. For example, a short range repulsion char- 
acterized by a length S will not significantly affect the 
free energy when average particle separations are larger 
than S. In these cases, the free energy will be governed 
essentially by entropy and our results should apply. For 
more confined systems, i.e. when particle separations are 
less than or comparable to 5, whether or not the free en- 
ergy is dominated by entropy or the potential will depend 
entirely on the details of the potential. We are currently 
investigating these scenarios. 



VI. DISCUSSION AND CONCLUSIONS 

We have introduced a model system composed of three 
hard disks confined to a circular corral. Though simple, 
this system exhibits caging and cage breaking behaviors, 
reminiscent of dense liquids, and allows us to exactly cal- 
culate a free energy landscape in terms of a single system 
coordinate. Respectively, caging and cage breaking can 
then be understood as the system becoming trapped in 
local energy minima and eventually being thermally ex- 
cited over an energy barrier. 

As the size of the system decreases, exploration of the 
configuration space is increasingly hindered by disk-wall 
and disk-disk interactions, though cage rearrangements 
do eventually occur for all e > 0. Arrhenius scaling de- 
scribes the transition time r between wells when the en- 
ergy barrier is large and the system is small. However, 
this scaling fails for smaller energy barriers, when diffu- 
sion becomes the limiting factor in relaxation. 

The coordinate(s) one uses to express the energy land- 
scape must capture the behavior one is studying [ToTJ2~3] . 
Those used in the text successfully describe motion from 
one cage to the other, however, there are any number 
of coordinates one could use that bear no relevance to 
caging. For example, a calculation of the collective ra- 
dius of gyration of the disks may yield the same result for 
configurations near cage breaking as for those far away. 
In the same respect, a particular value of h tells one very 
little about the actual configuration or spatial location 
of the disks. Furthermore, that caging can be described 
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VIII. APPENDIX 

Here we present a derivation of n(h) based on geomet- 
ric arguments. Because the distribution n{h) has even 
symmetry about h — 0, we directly treat only the cases 
where h > 0. All three disks are distinguishable and have 
the same diameter d. In simulations, we used d = 2, but 
for the sake of generality here we assume no specific disk 
size. The true radius of the corral is Rc = 3?* + e, but 
for simplicity, we define R = d + e as the radius of the 
corral accessible to the centers of the disks, as illustrated 
in Fig. [11] 

While Eqn. ([5| is general and exact, it can be further 
simplified for our system. Recalling that H(fi, r*2, 7%) in 
Eqn. ([5]) maps all rotationally equivalent states to the 
same h, we may consider only a single orientation of the 
system, such as the one in Fig. [TT] We constrain disks 1 
and 2 to always lie along the same horizontal line, L12, at 
a distance y above the bottom of the corral. For a given 
h, disk 3 may lie anywhere along a second horizontal line, 



8 




The lower limit y m i n corresponds to when disks 1 and 
2 are in contact at the bottom of the corral. Geometry 
yields that this occurs at 



Vd 



R 



d? 



(12) 



In this case, the upper limit 2/ max can correspond to 
L 12 = d or L 3 = 0. The correct value is the one that 
minimizes y max and therefore keeps all particle centers 
within the allowed region. Solving Eqns. Q and §W§ 
with these conditions yields, respectively, 



FIG. 11. The radius of the solid outer circle is Rc — R + r 
and is the physical boundary, but the centers of the confined 
disks can only exist a distance < 7? from the center, indicated 
by the dashed inner circle. 



min{2i? -h,R- 



(13) 



L 3 , at a distance y + h above the bottom of the corral. 
Thus, the number of states n{h) can be written as the 
product of integrals 



n(h) oc / dy I dx\ I dx 



dxz 



(8) 



and the problem becomes determining the appropriate 
limits of integration. For a circle centered at (0, R), we 
determine the lengths of chords L12 and L 3 to be 



L 12 = 2v/2yi?-y 2 
L 3 = 2^2(y + h)R~(y- 



hf 



(9) 
(10) 



There are four cases that must be considered to calculate 
n(h) correctly, each of which is geometrically distinct. 



B. 



Case B (h < d, X3 < x\ < X2) 




FIG. 12. Case B. 



A. Case A (h > d) 

The first of four cases that must be considered is shown 
in Fig. 11 This is the simplest case in that disk 3 never 
comes into contact with disks 1 or 2. Taking that the 
chord L12 begins at x = 0, the lower and upper bounds 
of X2 are, respectively x\ + d and Li2- The lower and 
upper bounds of x\ are and L12 — d. Disk 3 may exist 
anywhere along the chord L 3 independently of x\ and 
X2, so the lower and upper bounds for X3 are and L3. 
Therefore, 

Vmax L 3 L 12 -d L12 

riA oc J dy J dx% j dx\ J dx 2 

!/mi„ Xx+d 

'(£12 -d) 21 



dy L 3 



(11) 



When h < d, disk 3 may come into contact with either 
disk 1 or disk 2. In these cases, n(h) must be split into 
parts and computed in a slightly different manner. 

Three scenarios contribute to n{h): (x 3 < x\ < x 2 ) ; 
(xi < X3 < x 2 ) ; (x\ < x 2 < X3). The first and last of 
these scenarios are symmetric and contribute equally to 
n(h). The second scenario is presented in sections VIII C 
and lVIIIDI 

For now, we consider only cases where the x-coordinate 
of disk 3 is not between those of disks 1 and 2. Such a case 
is shown in Fig. |12| As stated, there are equally as many 
states where X3 < Xi < x 2 as there are for x\ < x 2 < x 3 . 
Therefore, we will calculate the multiplicity for only the 
states where £ 3 < x± < x 2l and finally multiply by two 
to obtain the contribution to n(h). Wc define w as the 
minimum horizontal separation between disks 1 and 3, 
given by 



w 



y/d 2 - h? 



(14) 



The minimum horizontal separation between disks 1 and 
2 remains d. 

We also define a quantity 5 = (L3 — L12) /2, as shown 



in Fig. 12 If we take the starting point of chord L12 
to be 0, then the lower and upper bounds on X3 are, 
respectively, [—SjL^ — d — w]. From Fig. 12 the limits 



of integration for X2 are found to be [x\ £12] and 
those of x\ are [X3 + w, L12 — d]. Including the factor of 
2 from the symmetric states, Eqn. (fsT) becomes 



i/max £l2— d— W L 12 —d L12 

us oc 2 J dy J dxa J dx\ J dx% 

Vxoln (Li2-L 3 )/2 x 3 +w x t +d 

3 



J/max L 12 — d Lx2 X 2 —'W 

fic oc J dy J dx\ J dx2 J dx^ 

J/min X-i+d Xt+W 

~(L 12 -d) 



dy (L 12 - dy 



G 



(17) 



The lower and upper limits of y are found in the same 
way as in case A. 



D. Case D (0 < h < V3d/2, Xi<x s < x 2 ) 



The limits of integration for y in this case are found 
from solving the equation 



L3 = w + d + 5 



(16) 



which after substitution of the various terms yields a 
quadratic equation in y. The limits of integration for 
y are the roots of this equation, such that the minimum 
root is y m i n and the maximum root is y max - 



Case C ( < h < d, xi < X3 < X2 




FIG. 13. Case C. 

Complementary to Case B, we now consider the case 
where (x± < X3 < X2) and y/3d/2 < h < d. As shown 
in Fig. 13 the minimum horizontal separation between 
disk 3 and one of the others is again w. Also, h is large 
enough that disks 1 and 2 may still contact, hence the 
minimum separation between disks 1 and 2 is d. 

Proceeding as before, the limits of integration for x±, 
X2 and X3 respectively are [0, L r2 — d], [rri + d, L12], and 
[x\ + w,X2 — w}. Eqn. (JsJ) becomes 




FIG. 14. Case D. 



The final case to be considered is also complementary 



to case B. In Fig. 14 disk 3 is again kept between 1 
and 2, only now h is such that disks 1 and 2 never come 
into contact. The horizontal separation between 1 and 
2 is always > 2w. Inspection again yields the limits for 
xi, X2, and X3 which are [0, L12 — 2wV, [X3 + w, L12], and 
[x\ + w,Li2 — w]. Therefore, Eqn. pi becomes 



J/max L 12 — 2W L 12 —W L 12 

Tin oc J dy J dx\ J dx3 J dx2 

Vrain X\+w x^+w 

?/max 

= J dy (L 12 - 2u;) 3 (18) 

l/min 

The lower limit y m - m comes from geometry, 



!-Jn 



R- \/R 2 - w 2 , 



(19) 



while, similar to case A, the upper limit y m ax is the min- 
imum value corresponding to either L 3 = or the maxi- 
mum root of L±2 = 2w, 
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Vn 



\{2R -h,R + \JR 2 - w 2 }. 



(20) 



Given the results of all four cases, the generalized 
method to compute n(h) is given by 



middle) and 9 w tt (bottom). At first glance, it would 
appear that n(0) = 2n(7r). However, one must consider 
how the number of states changes in the vicinity of and 
tt. Thus we are interested in 



n{h) 




h>d 

3d/2 <h<d 
0<h< y/3d/2 



(21) 



In general, the integrand for each case above is a cum- 
bersome function of h, thus we calculate n{h) using nu- 
merical integration. From Eqn. (21), the origin of the 



kinks in F(h) become more apparent as transitions from 
one solution regime to another. For d = 2, transition 
points occur at h = ±v3 and h = ±2, as described in 
the text. 



E. Behavior of Fb (h) as e — > 

To understand the growth of the energy barrier as e — > 
0, we only need the result of Case D. Setting h = gives 
w = d, and in the limit of small e, transitions occur along 
the diameter of the corral, which implies L 12 « 2c? + 2e. 



n(58) 

hm — - 

56->o n(TT — ot 



(25) 



In the top of Fig. 15 disks 1 and 2 are separated by 



an average distance s and disks 1 and 3 are separated by 
2s. The number of states where 9 ps is proportional to 
the product of the arc lengths subtended by 59 for disks 
2 and 3. Therefore, for the top configuration, n top oc 
(sS9)(2sS9) = 2s 2 59 2 . In the same way for the middle 
configuration, n m i<j = n top . Therefore, n(59) — 2n top oc 
As 2 S9 2 . For the bottom configuration where 9 « tt, the 
number of states is again proportional to the product of 
the arc lengths, but in this case the distance between 
disks 1 and 3 and disks 1 and 2 is s. Therefore, n(ir — 
59) = ribot oc s 2 59 2 . Thus, the ratio in Eqn. (25) is equal 
to 4. 



Eqn. ( 18 ) then becomes 



where 



n(h) oc / dy ( - ) e A oc e A Ay 



Ay = 2\/R 2 - d 2 = 2y / 2de + e 2 



(22) 



(23) 



In the limit e — >• 0, Ay grows as y/e and so n(h) cx e 7 / 2 . 
Therefore, the barrier height grows as 



7 

Fg oc — log n(h) oc — In e 



(24) 



as confirmed by the data shown in Fig. Qd) . 



2s 



0-0 



0-7T 




F. Barrier heights in 6 

Here we give an explanation as to why n(0) — 4n(7r). 
Shown in Fig. 15 are configurations where 9 ss (top and 



FIG. 15. Configurations where 6 w (top and middle) and 
6 rs tv (bottom). 
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